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We introduce a new discrete-time variational principle inspired by the quantum clock originally 
proposed by Feynman, and use it to write down quantum evolution as a ground state eigenvalue 
problem. This allows one to apply ground state quantum many-body theory to quantum dynamics, 
extending the reach of many highly developed tools from this fertile research area. We draw an 
explicit connection between previously known time-dependent variational principles and the new 
time embedded variational principle presented. Sample calculations are presented applying the idea 
to the electronic degrees of freedom in a Hydrogen molecule and the spin degrees of freedom of a 
model inorganic compound. We employ the configuration-interaction approximation on a basis in 
which time is a variable. Finally, we take advantage of the unique perspective of the new variational 
principle to examine the error of basis approximations in quantum dynamics. 



Introduction 

Feynman proposed a revolutionary solution to the 
problem of quantum simulation: use quantum comput- 
ers to simulate quantum systems. While this strategy 
is powerful and elegant, universal quantum computers 
may not be available for some time, and in fact, ac- 
curate quantum simulations may be required for their 
eventual construction. In this work, we will use the clock 
Hamiltonian originally introduced by Feynman [TJ [2] for 
the purposes of quantum computation, to re-write the 
quantum dynamics problem as a ground state eigenvalue 
problem. We then generalize this approach, and obtain a 
novel variational principle for the dynamics of a quantum 
system. Variational principles play a central role in the 
development and study of quantum dynamics [3HTU] , and 
the variational principle presented here extends the ar- 
senal of available tools by allowing one to directly apply 
efficient approximations from the ground state quantum 
many-body problem to study dynamics. 

Previous approaches for recasting quantum dynamics 
as a time-independent problem include Floquet theory 
for periodic potentials pTtfl3] and more generally the 
(£, t r ) formalism of Peskin and Moiseyev [H] . However, 
the approach proposed in this manuscript differs consid- 
erably from these previous approaches. We derive our 
result from a different variational principle, and in our 
embedding the dynamics of the problem are encoded di- 
rectly in its solution, as opposed to requiring the con- 
struction of another propagator. Moreover, the method 
of solution we describe will differ considerably from pre- 
vious approaches by employing techniques from ground 
state many-body theory to the time-dependent path as a 
whole. Considerable work has now been done in the mi- 
gration of knowledge from classical computing to quan- 
tum computing and quantum information [15 j. In this 
paper, we propose a novel use of an idea from quantum 
computation for the simulation of quantum dynamics. 

The successful application of quantum computers to 
quantum chemistry p~6j [17] has been demonstrated, in 



both theory [T8H2T] and experiment [22j[23]. The exact 
classical simulation of quantum dynamics is a challenging 
problem due to the exponentially large space in which 
the problem must be solved. To make progress using 
classical computers, accurate approximate schemes with 
controllable errors must be developed. Some ideas from 
quantum information have successfully led to new meth- 
ods for simulating many-body states on classical comput- 
ers [24j [25] . More generally, quantum information tools 
[T5] [26] have become widely applicable in a variety of set- 
tings, such as understanding the role of entanglement in 
many-body systems [27] . 

The paper is organized as follows. We will first review 
the Feynman clock: a mapping stemming from the field 
of quantum computation that can be employed for con- 
verting problems in quantum evolution into ground-state 
problems in a larger Hilbert space. We then generalize 
the Feynman clock into a time-embedded discrete varia- 
tional principle (TED VP) which offers additional insight 
to quantum time-dynamics in a way that is complemen- 
tary to existing differential variational principles [28]- 
[30] . The basic dynamics of an excited Hydrogen molecule 
are used to explain how dynamics is simulated with the 
clock Hamiltonian. We then apply the configuration in- 
teraction method [3TJ [32] from quantum chemistry to 
solve for approximate dynamics of a model spin system 
demonstrating convergence of accuracy of our proposed 
approach with level of the truncation. Finally we discuss 
metrics inspired by our approach that can be used to 
quantitatively understand the errors resulting from trun- 
cating the Hilbert space of many-body quantum dynam- 
ics. 



Physical dynamics as a sequence of quantum gates 

Consider a quantum system described by a time- 
dependent wavefunction \^{t)). The dynamics of this 
system are determined by a Hermit ian Hamiltonian H(t) 
according to the time-dependent Schrodinger equation in 



atomic units, 
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A formal solution to the above equation can be generally 
written: 

!*(*)> = T (e- if ?o dt ' H ^) |*(t )> = U(t,t ) |*(t )> 

(2) 

Where T is the well known time-ordering operator and 
U (£, to) is a unitary operator that evolves the system from 
a time to to a time t. These operators also obey a group 
composition property, such that if to < ti < ... <t n < t 
then 
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From the unitarity of these operators, it is of course also 
true that U{t n ,t n -\)^ — U(t n -i,t n ) where f indicates 
the adjoint. Thus far, we have treated time as a contin- 
uous variable. However, when one considers numerical 
calculations on a wavefunction, it is typically necessary 
to discretize time. 

We discretize time by keeping an ancillary quantum 
system, which can occupy states with integer indices 
ranging from |0) to \T — 1) where T is the number of 
discrete time steps under consideration. This quantum 
system has orthonormal states such that 



(4) 



This definition allows one to encode the entire evolution 
of a physical system by entangling the physical wave- 
function with this auxiliary quantum system represent- 
ing time, known as the "time register". We define this 
compound state to be the history state, given by 

i*> = ^X)i**)®i t > ( 5 ) 

where subscripts are used to emphasize when we are con- 
sidering a time-independent state of a system at time t. 
That is, we define |^) = \V(t)) \ t = ti . From these defini- 
tions, it is immediately clear from above that the wave- 
function at any time t can be recovered by projection 
with the time register, such that 
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Additionally, we discretize the action of our unitary op- 
erators, such that U(ti,to) = Uo and we embed this op- 
erator into the composite system-time Hilbert space as 
(Uo ® |1) (0|). While this notation is standard in the field 
of quantum information, it is less prevalent in chemical 
physics. For clarity, we offer an example algebraic ma- 
nipulation in the supplementary information. 

While the utility of this discretization has not yet been 
made apparent, we will now use this discretization to 
transform the quantum dynamics problem into a ground 
state eigenvalue problem. 
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FIG. 1. A schematic representation of the action of the clock 
Hamiltonian on the history state with three discrete times, 
and a Hilbert space of three states. Each block is a matrix 
with dimension of the physical system. 



Feynman's clock 

In the gate model [TSJ [33] of quantum computation, 
one begins with an initial quantum state |^o) state, ap- 
plies a sequence of unitary operators {Ui}, known as 
quantum gates. By making a measurement on the final 
state |^/), one determines the result of the computa- 
tion, or equivalently the result of applying the sequence 
of unitary operators {Ui}. The map from the sequence of 
unitary operators {Ui} in the gate model, to a Hamilto- 
nian that has the clock state as its lowest eigenvector is 
given by a construction called the Clock Hamiltonian [2] . 
Since its initial inception, much work has been done on 
the specific form of the clock, making it amenable to 
implementation on quantum computers [34]. However, 
for the purposes of our discussion that pertain to imple- 
mentation on a classical computer, the following simple 
construction suffices 



T-l 



H = Co + ( J ® I*) W " 0t ® |* + 1) (t\ 
1 t=o 

-U} <g> \t) (t + 1| + / \t + 1) (t + 1|) (7) 

where Co is a penalty term which can be used to specify 
the state of the physical system at any time. Typically, 
we use this to enforce the initial state, such that if the 
known state at time t = is given by |^o)- For this, we 
define: 



Co = (I-|*o)<*o|)®|0> <0| 



(8) 



One may verify by action of T-L on the history state 
defined above, |$), that the history state is an eigenvector 
of this operator, with eigenvalue 0. 
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A discrete-time variational principle 

We now introduce a new time-embedded discrete varia- 
tional principle (TED VP) and show how the above eigen- 
value problem emerges as the special case of choosing a 
linear variational space. Suppose that one knows the ex- 
act form of the evolution operators and wavefunctions at 
times t = 0,1. By the properties of unitary evolution it 
is clear that the following holds: 

2-(^i|/7 |^o)-(^o|^ t l^i)=0 (9) 

From this, we can see that if the exact construction 
of Ui is known for all i, but the wavefunctions are only 
approximately known(but still normalized), it is true that 

T-1 

J2 (2 - (* t+ i| u t |* t ) - <* t | u} |* t+1 >) > (10) 

£=0 

where equality holds in the case that the wavefunction 
becomes exact at each discrete time point. Reintroducing 
the ancillary time-register, we may equivalently say that 
all valid time evolutions embedded into the composite 
system-time space as \§ t ) \t) minimize the quantity 

S = J2(t'\(*t'\U\* t )\t) (11) 

t,t' 

where H (script font for operators denotes they act in 
the composite system-time space) is the operator given 

by 

T-1 

H=J2 I ®\ t ) (t\-U t ®\t + l) (t\ 

t=0 

-U} <8> \t) (t + 1| + 1 <g> \t + 1) (t + 1| (12) 

It is then clear from the usual ground state variational 
principle, that the expectation value of the operator 

s = J2(t'\(*t>\n\* t )\t) (13) 

t,t' 

is only minimized for exact evolutions of the physical 
system. This leads us immediately to a time-dependent 
variational principle for the discrete representation of a 
wavefunction given by: 

5S = 6Y,(t'\(*t>\H\* t )\t) =0 (14) 

t,t f 

It is interesting to note, that this is a true variational 
principle in the sense that an exact quantum evolution 
is found at a minimum, rather than a stationary point 
as in some variational principles [30]. This is related to 
the connection between this variational principle and the 
McLachlan variational principle that is explored in the 



next section. However, to the authors knowledge, this 
connection has never been explicitly made until now. 

To complete the connection to the clock mapping given 
above, we first note that this operator is Hermit ian by 
construction and choose a linear variational space that 
spans the entire physical domain. To constrain the solu- 
tion to have unit norm, we introduce the Lagrange mul- 
tiplier A and minimize the auxiliary functional given by 

C = £ (t'\ I H |* t ) \t)+C -X I 2 <*H*t> \t) - 1 

t,t' \t,t' 

(15) 

It is clear that this problem is equivalent to the exact 
eigenvalue problem on T~L with eigenvalue A. The op- 
timal trajectory is given by the ground state eigenvec- 
tor of the operator H. From this construction, we see 
that the clock mapping originally proposed by Feynman 
is easily recovered as the optimal variation of the TEDVP 
in a complete linear basis, under the constraint of unit 
norm. Note that the inclusion of Co as a penalty term 
to break the degeneracy of the ground state is only a 
convenient construction for the eigenvalue problem. In 
the general TEDVP, one need not include this penalty 
explicitly, as degenerate allowed paths are excluded, as 
in other time-dependent variational principles, by fixing 
the initial state. 

We note, as in the case of the time-independent varia- 
tional principle and differential formulations of the time- 
dependent variational principle, the most compact so- 
lutions may be derived from variational spaces that 
have non-linear parameterizations. Key examples of this 
in chemistry include Hartree-Fock, coupled cluster the- 
ory [35-37 , and multi-configurational time-dependent 
Hartree [38]. It is here that the generality of this new 
variational principle allows one to explore new solutions 
to the dynamics of the path without the restriction of 
writing the problem as an eigenvalue problem, as in the 
original clock construction of Feynman. 

Connection to previous time-dependent variational 
principles 

In the limit of an exact solution, it must be true that all 
valid time-dependent variational principles are satisfied. 
For that reason, it is important to draw the formal con- 
nection between our variational principle and previously 
known variational principles. 

Considering only two adjacent times t and t + 1, the 
operator T~L is given by: 

H=^(I®\t) (t\-U t ®\t + l) (t\ 

-U} <g> \t) (t + 1| + I <8> \t + 1) (t + 1|) (16) 

Now suppose that the separation of physical times be- 
tween discrete step t and t + 1, denoted dt is small, and 
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FIG. 2. The spin triangle within the vanadium compound, 
used as a model system for the TEDVP. Note that coor- 
dinating sodium ions and water molecules are not depicted 
here. The chemical formula of this compound is given 

by (CN 3 H6)4Na2[H 4 V^ V) 08(P04)4 (OCH 2 ) 3 CCH 2 OH 2 ] • 
14H 2 

the physical system has an associated Hamiltonian given 
by H, such that 

U t = e~ iHdt « I - iHdt - (17) 

By inserting this propagator into equation [l4j rewrit- 
ing the result in terms of and dropping terms of 
order dt 3 we have 

S - + dt)\ (I - iHdt) |*(t)> 

-(V(t)\(I + iHdt)\V(t + dt)) + 
{^(t)\H 2 dt 2 |*(t)) + (tf (t + dt)\V(t + dt))) = (18) 

After defining the difference operator such that 
d t \V(t)) = [\V(t + dt)) - \ V(t))] /dt, we can factorize the 
above expression into 

S (V(t)\ (idt - H) ] (idt - H) \V(t)) = (19) 

In the limit that dt —> 0, these difference operators be- 
come derivatives. Defining = idt\^(t)}, this is pre- 
cisely the McLachlan variational principle [30] . 

5||e-^|*(i))|| 2 = o (20) 

We then conclude that in the limit of infinitesimal phys- 
ical time for a single time step, the TEDVP is equivalent 
to the McLachlan variational principle. Under the rea- 
sonable assumptions that the variational spaces of the 
wavefunction and its time derivative are the same and 
that the parameters are complex analytic, then it is also 
equivalent to the Dirac-Frenkel and Lagrange variational 
principles [39 . Moreover, in an appendix, we provide a 
connection that allows other variational principles to be 
written as eigenvalue problems, and further discuss the 
merits of the integrated formalism used here. 

To conclude this section, we highlight one additional 
difference between practical uses of the TEDVP and 



other variational principles: The TEDVP is independent 
of the method used to construct the operator U t - In 
quantum information applications, this implies it is not 
required to know a set of generating Hamiltonians for 
quantum gates. Additionally, in numerical applications, 
one is not restricted by the choice of approximate prop- 
agator used. In cases where an analytic propagator is 
known for the chosen basis, it can be sampled explicitly. 

Structure of the TEDVP eigenvalue problem 

Suppose that the dimension of the physical system is 
given by N and the number of timesteps of interest is 
given by T. Assuming that the time register \T) is or- 
dered, the resulting eigenvalue problem is block tridiag- 
onal with dimension NT (See Fig. 1). 

This structure has been described at least once before 
in the context of ground-state quantum computation [40] , 
but to the authors knowledge, never in the context of 
conventional simulation of quantum systems. The block 
tridiagonal structure is particularly amenable to the use 
of parallel sparse eigensolvers. This opens up the possi- 
bility of exploiting the TEDVP to carry out simulations 
of time dynamics in parallel using modern computer ar- 
chitectures. Parallel-in- time physical dynamics [4TJ 02] 
has been seldom explored, but will become increasingly 
important as modern computer core speeds level out, and 
massive parallelism becomes the only way to accelerate 
computation. 

The solution of eigenvalue problems by sparse linear 
algebra methods is often limited by the gap between the 
lowest eigenvalue and the next largest eigenvalue. The 
spectral properties of the clock Hamiltonian have previ- 
ously been explored in ground-state quantum computa- 
tion [43] . Interestingly, it turns out that this Hamiltonian 
can be viewed as an analog to the tight-binding quantum 
particle in a box, and has an eigenvalue gap that scales as 
0(1/T 2 ), independent of both N and the physical time 
that the evolution of each operator U t represents. As the 
gap between the lowest and next lowest eigenvalue is a 
determining factor in the difficulty of applying iterative 
eigenvalue solvers, it is expected that the computational 
cost of this increase accordingly with the number of time 
steps taken, independent of the physical time each step 
represents. It is for this reason that we do not formu- 
late a continuous generalization of our variational princi- 
ple as in principle the number of time steps in this case 
tends to infinity. While the authors do not discount the 
formal value of a continuum construction, the computa- 
tional tractability of its practical implementations remain 
unclear. 

The eigenvectors of the clock Hamiltonian, which we 
will denote as time-orbit als, are interesting objects which 
merit study in their own right. The ground state time- 
orbital contains the physical dynamics of the system and 
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FIG. 3. The electron density of the first few clock eigenstates in an H2 molecule are shown at the 5 percent contour. The lowest 
state represents the density of the exact correlated dynamics, a slightly left-polarized Hi evolves for a half a period during 
the dynamics. Higher states show the localization into short spans of time in non-norm-conserving trajectories. Lobes pinched 
between the nuclei (black dots) in states 2 and 9 reflect significant double excitation character. With the applied penalty these 
states differ from one another by less than 3 milliHartrees. 



the excited time-orbit als contain orthogonal information. 
In order to examine the properties of these time-orbitals, 
we have diagonalized the clock representation of the full 
configuration interaction of the hydrogen molecule in a 
minimal basis (STO-3G [44]) beginning from an even su- 
perposition of the ground and first excited state. The 
electron densities of the resulting time orbitals are shown 
in Fig. 3. The lowest energy time-orbital, which provides 
the exact quantum dynamics given the initial condition 
conserves electron density. Interestingly the higher ly- 
ing states do not conserve probability over simulation 
time; instead probability localizes in time. The physical 
meaning of these states, and their time-reversal symme- 
try breaking is unclear and warrants further study. 

TEDVP for approximate dynamics 

There has been a recent rise in the interest of methods 
for simulating quantum spin dynamics in chemistry 
[45j [46] . To study the properties of the clock mapping 
when used to formulate approximate dynamics, we chose 
a simple model spin system inside an inorganic molecule 
[47] . Specifically, we examine the spin dynamics of the 
vanadium compound depicted in Fig. 2. 

By choosing the three unpaired electron spins to inter- 
act with one another by means of isotropic exchange as 
well as uniform static external magnetic field i?o, and a 
time-dependent transverse field B\, this system can be 
modeled with a spin Hamiltonian 

H = J a (5i • S 2 + Si • S 3 ) + J C S 2 • S 3 + 
fjLBoiSf + Si + Si) + fiB^Sf + S x 2 + Si) (21) 

Where Sf is the a Pauli operator on spin z, \i = g^, /i& 
is the Bohr magneton, g is the measured spectroscopic 



splitting factor. The couplings J a J c as well as g are 
fitted through experimental determinations of magnetic 
susceptibility [47 . The fact that they are not equal is re- 
flective of the isosceles geometry of the vanadium centers. 
The parameters of this model, are given by: g = 1.95, 
J a = 64.6 ± 0.5K, and J c = 6.9 ± IK. We will allow B 
to vary to study the properties of the clock mapping in 
the solution of approximate quantum dynamics. 

Dynamics in a linear basis 

The quantum chemistry community has decades of 
experience in developing and employing methods for 
obtaining approximate solutions of high-dimensional, 
ground- state eigenvector problems. By utilizing the con- 
nection we have made from dynamics to ground state 
problem, we will now borrow and apply the most concep- 
tually simple approximation from quantum chemistry: 
configuration interaction in the space of trajectories [48] , 
to our model problem to elucidate the properties of the 
clock mapping. 

For the uncorrelated reference, we take the entire path 
of a mean-field spin evolution that is governed by the 
time-dependent Hartree equations, and write it as: 

i*^>=E(n^i°>*)i*> ( 22 ) 
=E(ni°)i)i*) (23) 

= £lW> (24) 

t 

where 1 0) ^ is the reference spin-down state for spin z, |0)* 
is the reference spin-down state after rotation at time t, 
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FIG. 4. The convergence of the CI expansion in time is robust to strength of perturbation, and choice of initial state for the 
Vanadium complex. In these plots the error of the expectation value of S z at site 1 is plotted as a function of time, and coupling 
constant, each propagation with a different coupling constant is begun from a different random initial product state. 



\<pt) is the whole product system at time and U\ is 
determined from the mean field Hamiltonian. The equa- 
tions of motion that determine JJ\ are derived in a man- 
ner analogous to the time-dependent Hartree equations, 
and if one writes the wavefunction in the physical space 
as 

\m) = n u m i°>i = °(*) n i&> = °(*) i $ (*)> 

i i 

(25) 

Then the equations of motion are given by 

a(t) = a(0) (26) 
iU % = - E(t))Ui (27) 

Where if« is the mean field Hamiltonian for spin i 
formed by contracting the Hamiltonian over all other 
spins j ^ z, E(t) is the expectation value of the Hamil- 
tonian at time t, and / is the number of spins in the 
system. 

To generate configurations, we also introduce the 
transformed spin excitation operator S^, which is de- 
fined by 

S+lO^ll); (28) 
E/fS+t/f = S+ (29) 
^|0)* = |1)* (30) 

In analogy to traditional configuration interaction, we 
will define different levels of excitation (e.g. singles, dou- 
bles, ...) as the number of spin excitations at each time 
t that will be included in the basis in which the problem 
is diagonalized. For example, the basis for the configura- 
tion interaction singles(CIS) problem is defined as 

Beis = [Sf t \4> t ) \t) | j G [0,n],t e [0,T)} (31) 



Note that S^ t for j = n is simply denned to be the iden- 
tity operator on all spins so that the reference configura- 
tion is naturally included. Similarly, the basis for single 
and double excitations(CISD) is given by 

Bcxsv = {S+S+ \4> t ) \t) | j G [0,4*6 [0,j),te [0,T)} 

(32) 

Higher levels of excitation follow naturally, and it is clear 
that when one reaches a level of excitation equivalent to 
the number of spins, this method may be regarded as full 
configuration interaction, or exact diagonalization in the 
space of discretized paths. 

The choice of a time-dependent reference allows the ref- 
erence configuration to be nearly exact when Bq,B\ » 
J a ,Jc, independent of the nature of the time-dependent 
transverse field. This allows for the separation of the 
consequences of time-dependence from the effects of two 
versus one body spin interactions. 

After a choice of orthonormal basis, the dynamics of 
the physical system are given by the ground state eigen- 
vector of the projected eigenvalue problem 

H Bi \$)=E\$) (33) 

where we explicitly define the projection operator onto 
the basis Bi as P#. = Yl,\j)eBi b') 01 so ^ na ^ ^ ne P r °j ec ted 
operator is given by H Bi = {P&{HP&?) 

Using these constructions, we calculate the quantum 
dynamics of the sample system described above. For 
convenience, we rescale the Hamiltonian by the value of 
I//1B0. To add arbitrary non-trivial time dependence to 
the system and mimic the interaction of the system with 
a transverse pulse, we take B\ oc exp (— 1 2 /2) cos (rat). 
The magnitude of Bo was taken to be 200T in order to 
model perturbative two body interactions in this Hamil- 
tonian. To propagate the equations of motion and gener- 
ate the propagators for the clock operator we use the 
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FIG. 6. Trajectories for single-particle observables in differ- 
ent levels of basis set truncation starting from the reference 
configuration path (MF). As more levels of spin excitation 
are included with configuration interaction singles (CIS) and 
configuration interaction with singles and doubles(CISD) the 
trajectories converge to the exact result. 



FIG. 5. The total contribution to the ground state eigenvector 
from each level of excitation is plotted, where 0=MF, 1=CIS, 
2=CISD, and 3=FCI (Exact), and seen to decrease with exci- 
tation supporting the quality of the time-dependent reference 
state. The deviation of the wavefunction norm from unity re- 
sulting from projection is seen to decrease monotonically with 
level of excitation. 

enforced time-reversal symmetry exponential propaga- 
tor [49] given by 

U t = exp (~ijH(t + dt)j exp (-i^H(t)j (34) 

The dynamics of some physical observables are dis- 
played (Fig. 6) for the reference configuration, single 
excitations, double excitations, and full configuration in- 
teraction. The physical observables have been calcu- 
lated with normalization at each time step. It is seen 
(Fig. 4) that as in the case of ground state electronic 
structure the physical observables become more accurate 
both qualitatively and quantitatively as the configura- 
tion space expands, converging to the exact solution with 
full configuration interaction. Moreover in Fig. 5 we 
plot the contributions from the reference configuration, 
singles space, doubles space, and triples space and ob- 
serve rapidly diminishing contributions. This suggests 
that the time-dependent path reference used here pro- 
vides an good qualitative description of the system. As 
a result, perturbative and other dynamically correlated 
methods from quantum chemistry may also be amenable 
to the solution of this problem. 



In principle, approximate dynamics derived from this 
variational principle are not norm conserving, as is seen 
in Fig. 5, however this actually offers an important in- 
sight into a major source of error in quantum dynamics 
simulations of many-body systems, which is the trunca- 
tion of the basis set as described in the following section. 
Simulation based on the conventional variational princi- 
ples that propagate within an incomplete configuration 
space easily preserve norm; however the trajectories of 
probability which should have left the simulated space 
are necessarily in error. 



Norm Loss as a measure of truncation error 

Conservation of the norm of a wavefunction is often 
considered a critical property for approximate quantum 
dynamics, as it is a property of the exact propagator re- 
sulting from time-reversal symmetry. However, if norm 
conservation is enforced in a propagation which does not 
span the complete Hilbert-space, one must account for 
the components of the wavefunction that would have 
evolved into the space not included in the computation. 
It's not immediately clear how population density which 
attempts to leave the selected subspace should be folded 
back into the system without being able to simulate the 
exact dynamics. This problem is sometimes glossed over 
with the use of exponential propagators that are guar- 
anteed to produce norm-conserving dynamics on a pro- 
jected Hamiltonian. Some more sophisticated approaches 
adjust the configuration space in an attempt to mitigate 



8 



0.6 
0.4 
0.2 
0.0 



-0.2 



2 0.8 

X 



u 



CD 
O 



0.6 
0.4 
0.2 
0.0 



— S(MVP) - ■ SD(TEDVP) 
■ i S(TEDVP) — Ex(MVP) 
SD(MVP) - ■ Ex(TEDVP). 








4 6 
Time(a.u.) 



8 10 



FIG. 7. The dynamics of a spin observable at a short time 
step, dt = 0.01, are indistinguishable when generated with 
equations of motion determined by the time embedded dis- 
crete variational principle(TEDVP) (dashed lines) and the 
McLachlan variational principle(MVP) (solid lines). Results 
are shown for propagations restricted to the space of single 
excitations from the initial state(S), double excitations (SD), 
and the full space(Ex). The associated error metrics, A/i(t) 
for the TED VP (dashed) and J\f 2 (t) for the MVP (solid), also 
yield nearly identical results, displaying peaks correlated with 
qualitative deviations from the exact trajectory. 



the error [50] . 

This discrepancy is at the center of the differ- 
ence between the approximate dynamics derived from 
the discrete variational principle here and the ap- 
proximate dynamics derived from the McLachlan vari- 
ational principle such as the multi-configurational 
time-dependent Hartree method. Mathematically, 
this results from the non-commutativity of the 
exponential map and projection operator denned 
above. That is P Bi (t)T (V* / P B . {t) ^ 

T (exp (—i J dt' ' Pfi.( t t}H(t')Pfi.(jt\y) for a Hermitian op- 
erator H. In an approximate method derived from 
the McLachlan or any of the other differential time- 
dependent variational principles, the projection is per- 
formed on the Hamiltonian. As the projection of any 
Hermitian operator yields another Hermitian operator, 
the dynamics generated from the projection are guaran- 
teed to be unitary if a sufficiently accurate exponential 



propagator is used. Conversely, projection of a unitary 
operator, as prescribed by the TEDVP, does not always 
yield a unitary operator. Thus for an approximate ba- 
sis, one expects norm conservation to be violated, where 
the degree of violation is related to the severity of the 
configuration space truncation error. This leads us to 
define a metric of truncation error which we term the 
instantaneous norm loss. We define this as 



tfi{t) = l-\\P Bt U t P Bt \* t ) 



(35) 



where \i& t ) is always assumed to be normalized, which 
in practice means that a renormalization is used after 
each time step here. However, as we proved above, in 
the limit of a short time step, with dynamics generated 
by a Hamiltonian, the TEDVP must contain essentially 
the same content as the McLachlan variational principle. 
For this reason, we propose an additional metric which 
is given by 



N 2 (t) = \\(H(t) - P Bi H(t)P Bi ) \9 t ) \\ 2 dt 2 
= \\V(t)\* t )fdt 2 



(36) 
(37) 



Where H(t) is the physical Hamiltonian. This is moti- 
vated by appealing to the McLachlan variational princi- 
ple and substituting from the exact Schrodinger equation 
that idt = H\^ t ) where H is the full (non-projected) 
Hamiltonian. By defining V(t) = (H(t) - P Bi H(t)P B .) , 
we recognize this as a perturbation theory estimate of the 
error resulting from the configuration basis truncation at 
a given point in time. 

To examine the quality of these metrics and to better 
understand the consequences of the non-commutativity 
of the exponential map and projection, we return to the 
sample spin system. In this case, we choose a basis for 
the propagation space based entirely on the initial state, 
and do not allow it to change dynamically in time as be- 
fore. We perform simulations in the space of single exci- 
tations (S) from the initial state, double excitations (SD), 
and in the full Hilbert space (Ex). Dynamics from the 
TEDVP are generated by first building the exact prop- 
agator then projecting to the desired basis set while dy- 
namics from the McLachlan variational principle (MVP) 
are modeled by projecting the Hamiltonian into the tar- 
get basis set and exponentiating. A renormalization is 
used after each time step in the first case. Although 
one could perform the simulation with a timestep where 
timestep error is negligible, we remove this component of 
the calculation for this example by making the Hamilto- 
nian time-independent. This allows direct analysis of the 
effect of timestep on non-commutativity deviations. 

In Fig. 7 we show the dynamics for the two approx- 
imate truncation levels (S and SD) with both methods 
and their associated error metrics (A/i(t) and A/*2(t)). De- 
viations in the qualitative features of the observable oc- 
cur after even the first peak of the proposed metrics. In 
this particular simulation, the configuration interaction 
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with singles and doubles spans all but one state in the 
full Hilbert space. The lack of this one state results in 
the large qualitative errors present, associated with the 
first and subsequent peaks present in these error metrics. 
The impact of later peaks is more difficult to discern, due 
to error in the wavefunctions, which accumulates as the 
propagation proceeds. 

As predicted by the connection between the TEDVP 
and the McLachlan variational principle, while A/"i and 
A/*2 are not identical for each case, in the short time limit 
they yield extremely similar information, which is high- 
lighted in Fig. 7 displaying the resulting longer time 
dynamics for a time step of dt = 0.01. In Fig. 8, however 
we explore the effects of a significantly larger time step 
and begin to discern the result of the non-commutativity 
discussed here. Recalling that because the Hamiltonian 
is time-independent in this case, the propagator used is 
numerically exact in both instances, so this effect is not 
a result of what would be traditionally called time step 
error, resulting from intrinsic errors of an integrator. In- 
terestingly, it is observed that M\ (t) begins to decay to a 
nearly constant value. This occurs because the action of 
projection after exponentiation breaks the degeneracy of 
the spectrum of the unitary operator, resulting in eigen- 
values with norms different than 1. As a result, repeated 
action and renormalization of the operator is analogous 
to a power method for finding the eigenvector associated 
with the largest eigenvalue. This effect is exacerbated by 
taking long time steps. 



Conclusions 

In this manuscript, we introduce a new variational 
principle for time-dependent dynamics inspired by the 
Feynman clock originally employed for quantum com- 
putation. Unlike other previously-proposed variational 
principles, the proposed TEDVP approach involves the 
solution of an eigenvalue problem for the entire time 
propagation. This perspective allows for readily em- 
ploying many of the powerful truncation techniques from 
quantum chemistry and condensed-matter physics, that 
have been developed for the exact diagonalization prob- 
lem. The possibility of parallelize "in time" by building 
a large eigenvalue problem is important, as there are not 
many practical methods for achieving this in large-scale 
supercomputers. We introduce two novel error metrics 
for the TEDVP that allow one to characterize the basis 
approximations involved. We demonstrated the features 
of the method by simulating the dynamics of a Hydro- 
gen molecule and a molecular effective spin Hamiltonian. 
Further research directions include the use of other ap- 
proximate techniques for the time dynamics such as the 
use of perturbation theory [51 or coupled-cluster ap- 
proaches [52] . 
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FIG. 8. The dynamics of a spin observable at a large time 
step, dt = 0.5, shows significant differences resulting from 
the non-commutativity of projection and exponentiation in 
dynamics generated by the time embedded discrete varia- 
tional principle(TEDVP) (dashed lines) and the McLachlan 
variational principle(MVP) (solid lines). Results are shown 
for propagations restricted to the space of single excitations 
from the initial state(S), double excitations (SD), and the full 
space(Ex). The corresponding error metrics, Ni(t) for the 
TEDVP (dashed) and Af 2 (t) for the MVP (solid), differ con- 
siderably in this case. 
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Appendix: On the general construction of 
eigenvalue problems from dynamics 

We have presented one path for constructing eigen- 
value problems from quantum dynamics problems so far, 
however it is instructive to illuminate precisely which 
part of our procedure allowed this. To do this, we will 
slightly generalize the Floquet-type Hamiltonians and 
demonstrate that the time embedding was the crucial 
feature that allows one to recast a dynamics problem as 
an eigenvalue problem. This is in addition to the choice 
to work in an integrated framework, which we will show 
simply allows for a convenient choice of basis. 

Recall the definition of a Floquet-type Hamiltonian 
given by 



F(t) = H(t) - id t 



(38) 



If one considers a finite time evolution for a length of 
time T, this operator is Hermitian in the basis of Fourier 
functions states given by 



|* nj -(t)> = I*,-) |n) = !$,•> e 2 * int ' T 



(39) 



where |3>j) is a time-independent state of the physical 
system when considering the generalized inner-product 

(n'| \n) = [ T dt'e-^in't'/T^int'/T 

T Jo 

(40) 

This is the the generalized Hilbert space first introduced 
by Sambe [53] and generalized by Howland [54] . Because 
this operator is Hermitian in this basis, by noting the 
similarity to the Lagrange Variational principle 



SL = S 



f 

Jo 



dt <*(t)| F|*(t)) = 



(41) 



minimization of L on this linear basis of Fourier states 
yields a Hermitian eigenvalue problem. Thus the time 
evolution can be reconstructed by solving the full time- 
independent eigenvalue problem in this basis, or by con- 
structing a surrogate evolution operator as in the (£, t') 
formalism of Peskin and Moiseyev [14]. The use of 
Fourier basis states to express time dependence is natural 
given the form of the operator F. That is, matrix ele- 
ments of the derivative operator dt have a trivial analytic 
expression in this basis. This represents the same general 
idea we have been discussing here, which is to consider 
states in a system-time Hilbert space. However, as the so- 
lution of this variational problem will may yield a station- 
ary point rather than a true minimum [30] , ground state 
techniques are not appropriate for this particular oper- 
ator. Moreover, this operator is not in general Hermi- 
tian when considering arbitrary basis functions of time. 
For example arbitrary choices of plane waves not corre- 
sponding to the traditional Fourier basis will yield a non- 
Hermitian operator. The operator F' = (1/2) (F + F^), 



which has been used in the past for the construction of 
approximate dynamics [55], is Hermitian, however it still 
suffers from a pathology that the optimal solution repre- 
sents a stationary point rather than a minimum. How- 
ever, the operator 



G = F \ F= ( H (t) - id t )\H(t) - id t ) 



(42) 



is always Hermitian and positive semi-definite by con- 
struction. Thus one can expand the system-time Hilbert 
space in any linear basis of time, and the optimal path 
in that space will be the ground state eigenvector of the 
operator G utilizing the above generalized inner-product, 
assuming we have broken degeneracy by introducing the 
correct initial state. This can be viewed as an applica- 
tion of the McLachlan time dependent variational prin- 
ciple. From this, we see it is the expression systems in 
a system-time Hilbert space which allows for the eigen- 
value problem construction. Moreover, we note that this 
is not limited to the use of Fourier time basis states, and 
that many more expressions of time dependence may be 
utilized to construct an eigenvalue problem within this 
framework. 

One may ask then, what was the objective of work- 
ing in the an integrated formalism, defined using unitary 
operators rather than differential operators. To see this, 
consider evaluating the system at a point in time with 
the Fourier construction above. One has to expand the 
system in all time basis functions and evaluate them at a 
time t. When one usually considers time, however, they 
are thinking of a parametric construction where a number 
t simply labels a specific state of a system. Embedding 
into the system-time Hilbert space with this idea would 
be most naturally expressed as delta-functions. However 
matrix elements of d t on this basis can be difficult to con- 
struct. In the ancillary time system framework used in 
the TED VP, however, time is easily expressed as a dis- 
crete parametric variable. One might also consider the 
use of discretized derivatives in the operator G. How- 
ever, balancing the errors in numerical derivatives and 
the increasing difficultly of solving the problem with the 
number of simultaneously stored time steps can be prac- 
tically difficult. 
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